Reconstructing Patchy Reionization from the Cosmic Microwave Background 



o 
o 

(N 

X) 

(D 

0^ 



I 

o 



in 

oo 
O 



Cora Dvorkin^'^ and Kcndrick M. Smith'^ 
^ Kavli Institute for Cosmological Physics and Department of Physics, 
University of Chicago, Chicago IL 60637, U.S.A. 
^Enrico Fermi Institute, University of Chicago, Chicago IL 60637, U.S.A. 
^ Institute of Astronomy, University of Cambridge, CB3 OHA UK 
(Dated: February 9, 2009) 

We introduce a new statistical technique for extracting the inhomogeneous reionization signal 
from friture high-sensitivity measurements of the cosmic microwave background temperature and 
polarization fields. If reionization is inhomogeneous, then the optical depth to recombination will 
be a function r(n) of position on the sky. Anisotropics in r(n) alter the statistics of the observed 
CMB via several physical mechanisms: screening of the surface of last scattering, generation of new 
polarization via Thomson scattering from reionization bubbles, and the kinetic Sunyaev-Zel'dovich 
effect. We construct a quadratic estimator Tim for the modes of the r-field. This estimator separates 
the patchy reionization signal from the CMB in the form of a noisy map, which can be cross- 
correlated with other probes of reionization or used as a standalone probe. A future satellite 
experiment with sufficient sensitivity and resolution to measure the lensed B-modes on most of the 
sky can constrain key parameters of patchy reionization, such as the duration of the patchy epoch 
or the mean bubble radius, at the ~ 10% level. 



I. INTRODUCTION 

Upcoming generations of cosmic microwave back- 
ground (CMB) experiments will make precise measure- 
ments of secondary anisotropics on small scales (2000 < 
£ < 10000) in temperature and E-mode polarization and 
will also have sufficient sensitivity to measure secondary 
B-modes, for example the ones arising from gravitational 
lensing of the primary E-mode. 

Secondary anisotropics are caused by fluctuations in 
the distribution of matter after recombination and are 
generated by gravitational lensing, inverse Compton scat- 
tering of photons by the hot intracluster medium (the 
thermal Sunyaev-Zel'dovich effect), and the Doppler ef- 
fect produced by Thomson scattering of photons from 
radially moving electrons (the kinetic Sunyaev-Zel'dovich 
effect (kSZ)). 

Secondary anisotropics arc also generated during the 
inhomogeneous, or patchy, phase of reionization. Ana- 
lytic studies and numerical simulations suggest that the 
period of reionization was more complex than a sudden 
change in the ionization fraction [1-3]. Reionization is 
an inhomogeneous process and the inhomogeneities con- 
tribute to the small scale CMB temperature and polariza- 
tion anisotropics. In contrast to other secondaries, such 
as gravitational lensing [4-11] or the Sunyaev-Zel'dovich 
effect (SZ) [12, 13], the potential science returns from 
patchy reionization have not been extensively studied in 
the context of the CMB. 

There are currently almost no observational con- 
straints on the evolution of the ionization fraction during 
the epoch of reionization. The Gunn-Pcterson trough in 
high-rcdshift quasar spectra has been observed, showing 
that the transition from partial reionization to a fully 
ionized universe (xe ~ 1) took place at redshift z ~ 6 
[14-16]. Large-scale E-mode polarization measurements 
from five-year WMAP data show a total optical depth to 
recombination r = 0.087 ± 0.017. If reionization is as- 
sumed instantaneous, this would imply a transition red- 



shift Zrci = 11.0 ± 1.4 with 68% confldence [17]. Com- 
bining these two observations, we therefore have good 
indirect evidence for patchy reionization but no detailed 
information. 

Predictions from simulations tell us that in order to 
allow the contributions of the first stars [18], the begin- 
ning of reionization could go until z ^ 30. However, be- 
cause of a lack of experimental data, this prediction is still 
uncertain. The large-scale E-mode is sensitive to reion- 
ization history, but only contains information about the 
(spatial) average ionization fraction Xe{z) during reion- 
ization, not the size or morphology of the ionized regions. 
The WMAP data are not sensitive enough to place strong 
constraints on the reionization history beyond determi- 
nation of the total optical depth r, but future large-scale 
E-mode measurements from Planck or CMBpol can con- 
strain up to ~ 5 principal components or redshift bins in 

Xe{z) [19, 20]. 

The patchy reionization signal is small and, therefore, 
will be difficult to separate from the other secondaries us- 
ing the power spectrum alone. The patchy contribution 
to the temperature power spectrum is smaller than the 
sum of the contributions from lensing and low-redshift 
kSZ on all angular scales [21-26]. The thermal SZ signal 
from galaxy clusters leads to a larger signal than the one 
caused by the kSZ effect, but it can be separated from 
the other secondaries due to its non-blackbody frequency 
dependence. In polarization, currently favored models 
produce B-mode polarization power spectra with an am- 
plitude 0.01/xK) that is significantly lower than the 
B-modes coming from gravitational lensing {'^ 0.1/iK) 
[27-30]. 

In this paper, we will propose an estimator which iso- 
lates the patchy reionization signal in the CMB. If reion- 
ization is patchy, then the optical depth to last scatter- 
ing will be a 2D field r^m rather than a constant r. Our 
estimator will reconstruct each mode of this field, 
within statistical noise, using the small change in the 
CMB statistics which is induced by the mode. This con- 
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struction was inspired by the well-known lens reconstruc- 
tion estimator [4, 5, 7-9], which reconstructs each mode 
(j)im of the CMB lens potential in an analogous way. We 
will show that our new estimator Tim isolates the patchy 
signal, in the sense that its expectation value is simply 
the underlying field Tim, with no contribution from the 
Gaussian part of the CMB. 

In §11 and §111, we describe our semi-analytic mod- 
eling of reionization and its effect on the CMB. Wc will 
split the signal from patchy reionization into three effects: 
direction-dependent screening of the acoustic peaks, po- 
larization generated by Thomson scattering from reion- 
ization bubbles, and temperature anisotropy generated 
by radial motion of the bubbles (i.e. the kSZ effect). 

In §IV and §V we will construct the estimator T£,„, 
by analogy with the quadratic estimator for lens recon- 
struction. Wc will show that the first two effects from 
§111 (screening and Thomson scattering) have an alge- 
braic form where the quadratic estimator construction 
applies, but the third effect (kSZ) does not, because the 
kSZ anisotropy is proportional to the velocity field dur- 
ing reionization, and this field is not highly correlated 
with the primary CMB. Consequently, our estimator Tim 
is not very sensitive to the kSZ effect, but should be a 
near-optimal statistic for extracting the screening and 
Thomson signals. Although a complete treatment which 
formally extracts all the signal-to-noise is quite complex 
(§V B) we show that a simple estimator (§V C) contains 
all the information in practice. 

In §VI we demonstrate our estimator with Monte Carlo 
simulations. We make simplifying assumptions, such as 
simulating the lensed component of the _B-mode as if it 
were a Gaussian field. In §VII we present our forecasts. 
A future all-sky B-mode experiment can constrain key 
reionization parameters, such as the mean radius of the 
ionized bubbles or the duration of reionization, at the 
^ 10% level. It might be possible to improve this mea- 
surement making use of a delensing procedure or cross- 
correlating the r-map with large-scale structure. We con- 
clude in §VIII. 



We represent the H II regions by ionized spherical bub- 
bles of radius R with a log-normal distribution [1, 33], 
with a characteristic size given by R and width of the 
distribution given by trini?,: 
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We assume that the number density of bubbles fluctuates 
as a biased tracer of the large-scale structure with bubble 
bias b. For simplicity, we will assume that the bias is 
independent of the bubble radius. 

We take our fiducial model to be described by the fol- 
lowing set of parameters: 

{r. Ay, b, R, ainii} = {0.084, 19.0, 6.0, 5 Mpc, ln(2)} (4) 

In Fig. 1 we show the evolution of the ionization frac- 
tion in redshift (upper panel) and the radius distribution 
(lower panel) for our fiducial reionization model. 
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II. REIONIZATION MODEL 

Throughout this paper, we use a fiducial cosmology 
defined by the WMAP5+BA0-hSN parameters from [31]: 

{nbh^nmh\h,nA,T,A,7i,} ^ (i) 

{0.02265, 0.137, 0.701, 0.721, 0.084, 2.16 x 10'^ 0.96} 

For consistency with CAMB [32] , we model the ioniza- 
tion fraction in the following way: 
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where y{z) = (1 -I- z)'^' . yre and Aj^ are free parameters 
of the model. We parameterize the reionization model by 
taking Aj, and the total optical depth r to be the model 
parameters, and treat j/^e as a derived parameter. 



FIG. 1: Average ionization fraction Xe{z) (upper panel) and 
the bubble radius distribution P{R) (lower panel) in our fidu- 
cial reionization model. 

The line-of-sight integral for the optical depth during 
the inhomogeneous period can be written as 



r(n, z) = (TTnp.o 



rfz'(l 



H{z') 



-Xe(n, z'), (5) 



where H{z) is the Hubble parameter, is the Thomson 
scattering cross section, Up^o is the present number den- 
sity of protons, and Xe{n,z) is the ionization fraction in 
direction n at redshift z. 

Under the Limber approximation, the angular power 
spectrum of T(n) during the inhomogeneous period of 
reionization can be written as 
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PA..A..(x,fc = ^/x), (6) 
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where Axe is the ionization fraction fluctuation and x is 
the comoving distance from the observer along the line- 
of-siglit. 

The three-dimensional power spectrum Pax^Ax^ is a 
sum of two terms [27, 29, 34], corresponding to 1-bubble 
(16) and 2-bubble contributions (25) to the power spec- 
trum of the ionized hydrogen: 
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The 1-bubble contribution to the three-dimensional 
power spectrum is given by [29] 



-PAx.Ax.i^) 



x,il -Xe)[F{k) + G{k)], 



with the functions F{k) and G{k) defined by 
JdRP{R)[V{R)]^[W{kR)]^ 



m = 

G{k) = 



/ dRP{R)V{R) 
P(|k-k'|)F(fc') 
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where V{R) = AttR^ /i is the volume of the bubble, P{R) 
is the log-normal distribution of bubbles defined in Eq. 
(3), P{k) is the matter power spectrum, and W{kR) is 
the Fourier transform of a real-space tophat window func- 
tion with radius R, given by: 



W{kR) 
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[sin(fci?) - kRcos{kR)] 
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The first term in Eq. (8) is associated with the shot noise 
of the bubbles. The 2-bubble contribution is given by 

Pax^AxM = [(1 - X,) ln(l - x,)I{k) - x,f P{k), 

(12) 

where I{k) is defined by 
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and b is the bubble bias. 

Note that the 2-bubble contribution does not go to 
zero as Xe ^ 1, due to the x\P{k) term in Eq. (12). 
This term corresponds to fluctuations in the free elec- 
tron density which are due to matter fluctuations alone, 
independent of the bubble distribution. Since it comes 
from both the homogeneous and inhomogeneous periods 
of reionization, we will not consider the x\P{k) term to 
be part of the patchy signal in this paper, although we 
also find that this term is a small contribution to the 
total power spectrum (Fig. 2). 



III. REIONIZATION AND THE CMB 
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FIG. 2: Power spectrum CJ^ of the optical depth r(n) 
(solid/black curve) compared to the xiP{k) contribution (see 
Eq. (12)) to CJ'^ alone (dotted/red curve), in our fiducial 
reionization model. The x^P{k) term does not vanish as 
— > 1, and should not be thought of as part of the patchy 
reionization signal, but it is a small contribution to the power 
spectrum. 



1. Screening: temperature anisotropy and polariza- 
tion from the surface of last scattering are mul- 
tiplied by e~'^^"\ If reionization is inhomogeneous, 
the screening effect generates B-modes in polariza- 
tion because the optical depth T{n) to recombina- 
tion is a function of the line-of-sight direction n. 

2. Thomson scattering: new polarization is generated 
during the patchy epoch via scattering of the local 
CMB temperature quadrupole by ionized bubbles. 
This effect also generates B-modes [27-30]. 

3. kSZ: new temperature anisotropy is generated from 
the radial motion of reionization bubbles relative to 
the observer [22-26, 35]. 

In this section we will write unified expressions 
(Eqs. (19), (23) below) for the patchy signal in tempera- 
ture and polarization. The physics of each of these three 
effects will be encoded in "response fields" Ti and Ei, 
which will be studied in more detail in subsequent sec- 
tions. 

First, we introduce some notation. Throughout this 
section we will parameterize line-of-sight integrals by Xj 
the comoving distance from the observer. We denote the 
total optical depth to distance x along the line-of-sight 
n by: 



T(n,x) = axripfi 



^ dx ,^ . 
^a;e(n, x) 



(14) 



Inhomogeneous reionization produces three effects in 
the CMB: 



To study the patchy reionization signal in polarization, 
we write the observed polarization [Q ± iU) as a line-of- 
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sight integral: 
{Q±tU){n) 
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where aj„(n, x) are the temperature quadrupole mo- 
ments at each position in space and (±2^2™) are the 
spin ±2 spherical harmonics. Overdots denote deriva- 
tives with respect to conformal time. The integral for- 
mally runs from x = to x = 00, but it only receives 
nonzero contributions from the epochs of reionization 
(0 ^ ^ 20 in our fiducial model) and last scattering 
(700 < z < 1300). 

We write Xe(n, x) as its angular average plus a fluctu- 
ation term: 



Xe{n, x) = Xeix) + Aa;e(n, x) 



(16) 



Now let us expand {Q±iU){n) in powers of Axe, keeping 
zeroth and first order terms. A formal way of doing this 
is to write: 

{Q±iU){n) = {Q±tU)o{n) (17) 
+o-T?Vo / -^Axe(n,x)((9±iC/)i(n,x), 

where {Q ± iU)o is the polarization from recombination 
and homogeneous reionization, and the x-dependent re- 
sponse field {Q ± iU)i is defined by: 



(0±.a).(a,x) = /;<ix'2«|jfP 



(18) 
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= e--(^)5±,,(n,x)- / dx're--5p±i(n,x'), 
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where S/6t{x') in the first line denotes the fimctional 
derivative. (Note that the response field is an E-mode 

and so we will use the notation af^}'^^ when we write it 
in harmonic space.) 

The two terms in the second line correspond physi- 
cally to Thomson scattering and screening respectively. 
In the Thomson term, an Xe fluctuation generates new 
polarization which is proportional to the source term S^^^ 
evaluated at the same point along the line-of-sight. In the 
screening term, an Xg fluctuation changes the amplitude 
of all polarization generated at earlier times, i.e. from 
earlier in the patchy epoch or from recombination. 

In this paper, we will find it convenient to represent 
Eq. (17) in binned form. We introduce a set of TV red- 
shift bins which cover the patchy epoch. We denote the 
endpoints of the a-th bin by [z,";^, z^ax] the central 
redshift by z". We denote the values of x which corre- 
spond to these values of z by Xmin' X": Xmax- (Through- 
out this section, we have found it convenient to take x 
as the time coordinate, but in subsequent sections it will 
be more intuitive to use z, so we define our binned rep- 
resentation using redshift bins instead of x bins.) 



We then write: 
(Q ± iU){n) = (0 ± tU)o{n) + ^ Ar"(n)(Q ± iU)f{n), 

a 

(19) 



where: 



Ar"(n) = arnp,o ^Axe(n,x) (20) 



(g±zC/)?(n) = e-^(>^°)5±i(n,x") 

'pol 



dx'fe-^S^Jn,x') (21) 



We have assumed here that AT"(ri) <c 1. 

Eqs. (19) and (21) are the representation for the 
patchy polarization signal that we will use throughout 
this paper. In this form, we think of patchy reionization 
as being represented by a set of 2D fields At" (n) which 
correspond to the optical depth anisotropy in redshift 
bins. The CMB polarization signal is obtained by multi- 
plying each of these fields by a "response field" {QztiU)f, 
which is defined (Eq. (21)) in a way which does not de- 
pend on the bubble realization. The response field is a 
pure E-mode; we will use the notation af^^"'^ when we 
write it in harmonic space. 

The temperature case is analogous but the source term 
is more complicated. We give an outline of the calcula- 
tion here; the details can be found in App. B. We can 
still write T{n) as a linc-of-siglit integral: 



/•OO 

T(ri)= / dxST{n,x,r{n,x)), 



(22) 



where St is a function of local quantities at (n, x) and 
the local optical depth r(n, x). (The full expression for 
St is given in Eq. (B6).) Starting from this line-of-sight 
integral, wc obtain a binned expression: 

r(n) = To(n)+^Ar"(n)Tf(n) (23) 

a 

, SSt 



Tf(a) 



(24) 



In App. B, we evaluate the functional derivative and 
compute the explicit form of the response field Ti (see 
Eq. (Bll)). The main contributions to Ti are a Doppler 
term in which a r fiuctuation generates new temperature 
anisotropy at the same redshift, and a screening term 
in which a r fluctuation changes the amplitude of all 
anisotropies generated at earlier redshifts. 

For now we simply note that the patchy signal in tem- 
perature can be written in a form (Eq. (23)) which is 
analogous to polarization (Eq. (19)): the optical depth 
(At") in each redshift bin multiplies a response field Ti 
which is independent of the realization of patchy reioniza- 
tion. We will use this representation for the temperature 
signal throughout the paper. 

Since we have split the Ar-ficld into redshift bins, we 
will need to split the power spectrum CJ'^ into bins as 
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FIG. 3: Comparison between the B-mode lensing power spec- 
trum (solid curve) and the B-mode power spectrum from 
patchy reionization (dotted curve) in our fiducial model. 



well. In the Limber approximation the power spectrum 
is diagonal in redshift bins, and given by restricting the 
integral in Eq. (6) to the relevant redshift range: 



C] 



p,0 



(25) 

In polarization, inhomogeneous reionization generates a 
B-mode, in contrast to the homogeneous case. Under the 
Limber approximation, valid for ^ ^ 10 [29, 30], the B- 
mode power spectrum is related to the power spectrum 
of the ionized hydrogen, PAx^Ax^ik) as 
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xe-2-(>^)g?^,(x)PA..A..(x, k = tlx). (26) 
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where Pax^Ax^ is given by Eqs. (7), (8) and (12), and Xii 
Xf are comoving distances to the beginning and end of 
reionization respectively. Qrmsix) is the rms temperature 
quadrupole, which we treat as constant during the epoch 
of patchy reionization, with a value of Qrms = 22/iK in 
the fiducial model. In Fig. 3 we show the B-mode power 
spectrum from patchy reionization (Eq. (26)) with the 
lensing B-mode shown for comparison. 

Throughout this paper, CMB noise power spectra will 
be given by 



EE 
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where Ap is the detector noise in /zK-steradians and 
^'fwhm is the FWHM of the beam in steradians. 



IV. QUADRATIC ESTIMATOR FOR PATCHY 
REIONIZATION: TOY EXAMPLE 



In the previous subsection, we saw that the CMB 
anisotropy with the patchy reionization contribution in- 
cluded can be written in the following form: 



Tin) = Toin) + J2Ar"{n)T,"in) 



(28) 



(Q ± iU)in) = (Q ± iU)o{n) + ^ AT"(n)(Q ± iC/)?(n) 

a 

(29) 

The fields {At)°',T^,{Q ± iU)'^ are not directly ob- 
servable, but their presence alters the statistics of the 
CMB. The leading terms Tq, {Q±iU)r, in Eqs. (28), (29) 
are Gaussian to a good approximation (ignoring non- 
Gaussian contributions from CMB lensing which will be 
discussed later) but the patchy terms are not Gaussian 
fields. 

How can we construct an estimator which can detect 
the patchy terms statistically? In this section and the 
next, we will propose an answer to this question. We will 
construct an estimator Tim which reconstructs the modes 
of the T field from the CMB. This estimator will separate 
the patchy and Gaussian contributions, in the sense that 
the Gaussian part of the CMB does not contribute to the 
expectation value of the estimator (but does act as noise 
in the reconstruction). 

In this section, we will consider a simplified case to 
build intuition. We consider polarization only and make 
the "constant quadrupole" approximation: we assume 
that {Q ± iU)i = (Q ± iU)^ for all pairs of redshift bins 
a, (3 during the epoch of patchy reionization. Physically, 
this approximation means that the source term for CMB 
polarization. 
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i^2m(n))aL(n,z) 



(30) 



depends on the line-of-sight direction n but is inde- 
pendent of the redshift z. In reality, the temperature 
quadrupole a|^(n, z) is not 100% correlated between red- 
shifts, but decorrelates on the scale Az k, 5 [30]. We will 
see how to include this decorrelation in our formalism in 
the next section. 

In the constant quadrupole approximation, Eq. (29) 
becomes: 

(Q±^c/)(n) = (g±ic/)o(a) + AT(fi)(g±zc/)i(a), (31) 

where Ar(ri) = ^^Ar"(n) is the inhomogeneous part 
of the total optical depth to recombination. There is a 
formal analogy between this expression and CMB lensing 
[6, 36, 37]. The lensed CMB polarization is given by: 

(Q±iC/)le„.cd = (Q±iC/)u„l+(V0)aV''(Q±ii7)„nl, (32) 

where neither the unlensed CMB {Q±iU)-an\ nor the lens 
potential (j) are directly observable. Nevertheless, there 
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is a quadratic estimator which "reconstructs" the 
lens potential mode by mode, using only the observed 
polarization {Q ± iC/)ienscd [4, 5, 7-9]. 

On a technical level, the lens reconstruction estima- 
tor (j)ern is possible because the two-point function of the 
lensed CMB contains a term which is linear in (h: 
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(33) 



where the lensing coupling F^^J^^'' is given in App. A 
(Eq. (A9)). (The expectation value (•) in Eq. (33) is 
taken over realizations of the CMB in a fixed realization 
of the lens potential (p.) 

Analogously, in the constant quadrupole approxima- 
tion of this section, patchy reionization induces a two- 
point function which is linear in the field At: 



£m ^ ' 

where the patchy coupling E^^Jj'' is given by: 
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Starting from this two-point function, one can construct 
a minimum- variance quadratic estimator r^™ for the field 
Ar; 



(36) 
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where 
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The derivation of this estimator is given in App. A. 
The estimator Tim is unbiased, in the sense that: 

(nm) = (Ar)fo„, 



(37) 



(38) 



and N^'^ is the noise power spectrum of the reconstruc- 
tion. 



{'''em^C m') noise — NJ"^Su'S„ 



(39) 



where (•) noise denotes an expectation value taken over 
CMB realizations containing no patchy signal (i.e. Gaus- 
sian realizations). The expectation value on the LHS 
of Eq. (37) is what we mean by "reconstruction noise" 
power in the context of r reconstruction; it is the power 



spectrum of the reconstruction due to statistical fluctua- 
tions alone, in the absence of any patchy signal. 

The construction of r^m in Eq. (36) is formally identical 
(with a different F coupling) to the lens reconstruction es- 
timator (pim] details are given in App. A. However, there 
is one qualitative difference between the patchy quadratic 
estimator Tim and the lensing estimator (j)i„i that wc 
would like to emphasize. The polarization signal from 
lensing is proportional to the unlensed CMB (Eq. (32)) 
whereas the patchy signal is proportional to a new field 
El (Eq. (31)). The quadratic estimator T£m depends 
on the fields Eq and Ei being correlated. If the cross 
power spectrum C^°^^ were zero, then the quadratic re- 
construction would not be possible (the reconstruction 
noise in Eq. (37) would be formally infinite). 

In Fig. 4 we show the cross power spectrum C^"^^, 
with the auto power spectrum Cf"^" shown for compar- 
ison. The cross power spectrum was calculated using a 
modified version of the CAMB code and the details of 
our calculation of Cf°^^ are given in App. B. 
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FIG. 4: Auto power spectrum Cf°^° and cross power spec- 
trum Cf^-^^^"' for varying redshift 9 < z < 21. Here, Eq is the 
E-mode without patchy reionization and Ei{z) is the response 
field for T fluctuations at redshift z, defined in Eq. (18). On 
large scales, Cf°^^ is positive and dominated by the Thomson 
contribution to the response field. On small scales, Cf°^^ is 
dominated by the screening contribution to the response field, 
and Cf"^i ^ -Cf°'^'>. 



Note that Cf°^^ is positive at large scales due to 
Thomson scattering, while at small scales is negative due 
to the e""^ screening. To understand this effect one can 
think of the response field Ei as the variation of the E 
mode with t fluctuations. Then, the addition of r fluctu- 
ations generates more power in the reionization peak, but 
the e"'^ screening becomes larger too, so the amplitude 
of the acoustic peaks is suppressed. 

In Fig. 5 we show a toy example of our estimator r^™ 
applied to a simulated polarization map (which was done 
using the HEALPIX package [38]), with the signal-to- 
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noise artificially increased by omitting the lensed B-mode 
from the simulation. (Note that the inhomogeneous r 
fluctuations will still generate a i3-mode). It is seen that 
the reconstruction Tim recovers the input field At from 
the observed CMB, within some statistical noise associ- 
ated with the reconstruction. 

The construction of the quadratic estimator Tim given 
in this section illustrates the qualitative features of our 
method, but depends on an approximation that we would 
like to relax: the line-of-sight CMB quadrupole is con- 
stant throughout the patchy epoch. Additionally, our es- 
timator only uses the EB two-point function (Eq. (34)). 
In principle extra signal-to-noise can be obtained by con- 
sidering all cross-correlations of {T,E,B}. For this rea- 
son, we will next consider a more general construction 
(§V) before presenting signal-to-noise forecasts in §VII. 



V. QUADRATIC ESTIMATOR FOR PATCHY 
REIONIZATION: COMPLETE TREATMENT 

A. Preliminaries 

In the previous section, we constructed a quadratic es- 
timator T£m for the optical depth anisotropy from patchy 
reionization, under the simplifying assumption of con- 
stant quadrupoles along the line-of-sight. In this section, 
we will consider an extended version of this construction 
which does not rely on this approximation, and which 
uses TT, TE, EE, and TB cross-correlations in addi- 
tion to EB. This will introduce significant additional 
complexity, but in the end we will show that a simple es- 
timator which is very similar to the previous one contains 
almost all of the signal-to-noise in practice. 

As shown in §111, if we divide the epoch of patchy reion- 
ization into redshift bins a, then the observed CMB can 
be written: 

T{n) = To(S) + 5]Ar"(S)Tf(S) (40) 

a 

iQ±iU){n) = (Q±zC/)o(a) + ^Ar"(S)(Q±z[/)?(fi) 

a 

All calculations in this section have been done using bins 
with Az = 0.6. 

From Eq. (40), it follows that the two-point function 
of the CMB, to first order in At, can be written: 



-m2 



(41) 
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where we have defined 
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where the F couplings are given by: 



(47) 

Note that the redshift dependence of the F couplings en- 
ters through the cross power spectra cf°^^'^"'\ since the 
response field Yi is redshift-dependent. In Fig. 6 we show 
each of these cross power spectra as a function of £ and z. 
The cross power spectra were calculated using a modified 
version of CAMB (App. B). 

To get some intuition into these cross spectra, first note 
that for high ^, the screening effect dominates and the re- 
sponse fields satisfy Ti « —Tq, Ei —Eq. Each cross 
spectrum Cf°^^ is approximately equal to {—Cf°^°). At 
low the Doppler contribution to Ti and the Thomson 
contribution to Ei become important, and one can see 
additional features in the power spectra. The largest ones 
are the low-£ bump in Cj°^^ (Fig. 6, lower panel), which 
arises from correlation between the Doppler contribution 
to To and the Thomson contribution to Ei , and the bump 
in Cf°^^ (seen previously in Fig. 4) arising from correla- 
tion between the reionization E-mode and the Thomson 
contribution to Ei. 



B. Principal component analysis 

In this subsection, we will construct quadratic estima- 
tors which formally extract all the signal-to-noise, given 
the TEB two-point function in Eq. (41). This will gener- 
alize the single-field EB estimator from §IV. In the next 
subsection, we will show how to simplify the construction 
by constructing a simple estimator which extracts nearly 
all the signal-to-noise in practice. 

In App. A, we show how to construct quadratic esti- 
mators from the TEB two-point function (Eq. (41)) using 
a principal component construction. The key results can 
be summarized as follows. 

There are N "eigcnmode" quadratic estimators 
defined by: 
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FIG. 5: An exaggerated example of the quadratic estimator Tern defined in §IV assuming the constant quadrupole approximation. 
Top row (left to right): the primary E-mode _Eo, the response field Ei and the Ar-field, on a 100 deg^ patch of sky. For visual 
purposes, the Ar-field has been multiplied by a Gaussian window function, and we have artificially increased the signal-to-noise 
of the reconstruction by omitting the lensed B-mode and assuming cosmic variance limited measurements to ^max ~ 2000. Note 
that El « —Eo on small scales, while the -Bi-field has power added at large scales. Bottom row (left to right): the _E-mode and 
B-mode components of the total observed polarization (Eq. (31)), and the quadratic reconstruction rim- Note that _B-modes 
appear where there are r fluctuations. In this figure, the units of the E-mode and B-mode polarization fields are in fiK and 
At is dimensionless. 
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^^^W 2(2^+1) 



(49) 



)ei )i2 ^ iii2i ^ tit 



-l^XX' ( r^^l-^YY' j.XY(t^)*j.X'Y' (rp) 

XX'YY'tit2 

where X,Y e {T, E, B} and 



'TT I atTT (jTE 



CP +N, 







TE (jEE ^ ^EE 



cp^ + nP^ 



(50) 

The redshift weighting w[^\z) which appears in the 



estimator is computed using a principal component anal- 
ysis described in App. A 2. In principle, the redshift 
weighting is ^-dependent, but in practice the weighting 
is independent of £ to a good approximation (Fig. 7). 
Note that the first eigenmode peaks at the redshift where 
Xe ~ 0.5 in our fiducial model, where we expect to get 
the largest contribution to the r-power spectrum. 

The noise covariance of the quadratic estimators (gen- 
eralizing Eq. (39) from the previous section) is given by: 

{£^ll f£9^/)noiso = N^^^'^ SijSil'Smm' , (51) 

i.e. the quantity A^J^^*'' defined in Eq. (49) is the re- 



9 




FIG. 6: Auto and cross power spectra Cf""^", Cfo^°, Cj°'^\ 
Cj°^^ and Cf""^^ which appear in the F couphngs (Eqs. (42)- 
(46)) and in the quadratic estimator. Notation follows §111: 
To,Eo are the temperature and E-mode polarization fields 
with patchy contribution omitted, and T\,E\ are the response 
fields to r fluctuations at redshift z. 



FIG. 8: Cumulative signal-to-noise versus maximum multi- 
pole I in the r reconstruction, for the eigenmode quadratic 
estimator (Eq. (48)) split into contributions from EB, TB, 
TT, EE, and TE cross-correlations. In both panels, we take 
= 3000 and Ap = V2At. In aU cases, the EB pair has 
the dominant contribution to the signal-to-noise. 
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FIG. 7: Redshift weighting «j^''(z) for the first and sec- 
ond quadratic estimators (i = 1,2) defined in Eq. (48). 
The redshift weightings have been normalized so that 

E„«'f'(^") = 1, and j:ci^?\nf = 1/16. It is seen 
that the weightings are nearly independent of £. 



construction noise power spectrum of the z-th quadratic 
estimator and the quadratic estimators are uncorrelated 
with each other. 

The principal component analysis that is used to con- 
struct the estimators S^^^^ also guarantees that the esti- 
mators remain uncorrelated even in the presence of the 
At signal. We can therefore compute the total signal-to- 
noise of the detection of patchy reionization by summing 
the signal-to-noise per ^-mode for the i-th estimator (de- 



noted by A^*'' and defined precisely in Eq. (A19)) over all 
values of i and £: 

S/N = 



^^(2£+l)(A«f 



1/2 



(52) 



where /sky refers to the observed fraction of the sky. 

In Fig. 8, we show a signal-to-noise forecast for two 
experiments. For the first experiment (Ap = 2 fiK- 
arcmin, ^fwhm = 10 arcmin), we find fsi^y{S/N)^ = 21. 
This experiment is mainly Hmited by the beam size; at 
this resolution patchy reionization is detectable only by 
a large- /sky survey such as a future CMB polarization 
satellite. (Some ideas for improving the S/N are given 
in §VIII.) For the second experiment (Ap = 0.5 fiK- 
arcmin, 6'fwhm = 1 arcmin), we find /sky (S'/iV)^ — 217. 

In Fig. 8 we also show the contribution to the total 
signal-to-noise from each of the 5 pairs of fields TT, TE, 
EE, TB, EB in the quadratic estimator. It is seen that 
almost all the signal-to-noise comes from the EB estima- 
tor. Also, note that the scale at which the signal-to-noise 
curves saturate {£ w 400) is given by the scale at which 
the signal peaks. 

Another way of splitting up the total signal-to-noise is 
to consider the contribution from each of the quadratic 
estimators. When this is done, it is seen (Fig. 9) that 

the first estimator £j,J^ contains almost all the signal- 
to-noise. Intuitively, this means that our reconstruction 
of patchy reionization has "the degrees of freedom of a 
single 2D field". (The opposite limiting situation, i.e. 
many quadratic estimators with roughly the same signal- 
to-noise but orthogonal redshift weightings, would corre- 
spond to reconstructing the degrees of freedom of a 3D 
field.) 
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C. A simple estimator which contains all the S/N 
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In the previous subsection, we constructed a complete 
set of eigenmode estimators which formally extract all 
the signal-to-noise given the TEB two-point function in 
Eq. (41). However, the conclusion of this analysis was 
that almost all the signal-to-noise was contained in the 
first principal component S^n- Intuitively, this means 
that all the redshift bins are so highly correlated that 
the choice of redshift weighting is not important. Fur- 
thermore, we found that all the signal-to-noise was con- 
tained in the EB estimator. 

These observations motivate the following construc- 
tion. We fix a redshift bin fi and construct an estimator 
by simply repeating the "toy" example from §IV, 



' f.m 



using the T couplings corresponding to bin /i: 



FIG. 9: Signal-to-noise per mode A^*' for the first three (i — 
1, 2, 3) eigenmode quadratic estimators f^'^, keeping only the 
EB contribution. It is seen that the first eigenmode contains 
nearly all the signal-to-noise. 
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(53) 



An interesting result emerges when we extend the anal- 
ysis to smaller scales (i.e. -^max > 3000). At those 
scales, the TT estimator acquires large signal-to-noisc. 
Fig. 10 shows the signal-to-noise coming from the TT 
pair as a function of ^max- In principle, at smaller scales 
(^max ^ 5000), the TT estimator is more sensitive than 
EB. However, it is unclear whether other secondaries such 
as point sources, kinetic SZ and thermal SZ effects will 
permit a clean measurement of patchy reionization on 
these scales. We defer a complete treatment of this issue 
to future work. 
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where: 



N 



2£ 



+ 1 ^ 



.SS(v)|2 



EE 



BB 



(54) 

We will now compute the signal-to-noise of this estimator 
in a way which does not rely on the constant quadrupole 
approximation from §IV, and show that the signal-to- 
noise is nearly the same as the optimal estimator from 
the previous subsection. 

Without the constant quadrupole approximation, the 
expectation value of the reconstruction has a nontrivial 
redshift dependence: 



(55) 



where the redshift response is given by 



n: 



r(M) 



p£;S(v)-p£;S(r„)* 

^£l£2f (1(21 



(56) 



The total signal-to-noise is given by 



S/N 



/sky 

2 



^(2^+1) % 




1/2 



(57) 



FIG. 10: Total signal-to-noise for the TT quadratic estimator, 
assuming cosmic variance limited temperature measurements 
up to some maximum multipole £max, and treating the lensing 
contribution to the temperature power spectrum as a source 
of noise. For £max ^ 4000, the signal-to-noise rises sharply 
but in practice the measurement is likely to be contaminated 
by non- Gaussian secondary anisotropies. 



where 



(58) 



If we consider a survey with Ap = 0.1 /xK-arcmin, 
?FWHM = 1 arcmin, £niax — 3000, and take the redshift 
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FIG. 11: Redshift response (z) for the simple estimator, 
defined in Eq. (56), sfiowing a near-unit response {Rg'^\z) ~ 
1) for all redshifts. This plot was made for a redshift bin fi 
that corresponds to z = 10.2. 



The forecasts in §VII will show that we expect a < lOa 
detection for futuristic sensitivity levels. The uncertainty 
ct(t) would have to be > 0.01 (i.e. a fractional uncer- 
tainty of 10%) in order to bias this measurement. The 
forecasted uncertainty from Planck is (t(t) w 0.0045 [39], 
so parameter uncertainty is unlikely to be an issue. 

In conclusion, the principal component construction 
from the previous subsection should be regarded as a 
proof that the simple estimator presented here extracts 
all the information. For practical analysis, the principal 
component analysis is not necessary and one can use the 
simple estimator constructed in Eq. (53). For forecasting, 
one can compute a noise power spectrum for the simple 
estimator using Eq. (54), and the signal power spectrum 
is given to an excellent approximation by Eq. (59). Un- 
certainty in cosmological parameters is unlikely to be a 
complicating factor. 



bin fi which corresponds to 2 = 10.2, we find the following 
numerical results. (The qualitative conclusions do not 
depend on these choices.) The signal-to-noise for the 
simple estimator (Eq. (57)) is nearly the same (^99.9%) 
as the principal component estimator (Eq. (52)). The 

redshift response Ri{z) is wl, i.e. the estimator t^^^ 
simply reconstructs the sum of the r fluctuations over all 
redshift bins (Fig. (11)). Intuitively, this is because the 
redshift bins are so highly correlated that the estimator 
which is optimized for t fluctuations in bin fi will pick up 
the fluctuations from every bin. The near-unit redshift 
response also implies that 

a. 

Numerically, we find that tliis agreement is better than 

1%. 

Let us conclude this section by discussing the issue of 
model uncertainty. Our construction of Tim depends on 
knowing the cross spectrum cf°^''^^\ which depends on 
cosmological parameters. If we have imprecise knowledge 
of the cosmology, does this affect the estimator in a signif- 
icant way (for example, by changing the normalization)? 

We expect that the most important parameter depen- 
dence will be the r dependence of the cross spectrum 
Cf°^^. (The dependence of the reionization E-mode on 
geometric distances, for example, should be much weaker 
than the r dependence.) At low the cross spectrum is 
roughly proportional to r, because the largest contribu- 
tion is from cross-correlating the patchy Thomson sig- 
nal with the primary reionization E-mode. At high 
the cross spectrum is roughly independent of r, provided 
that T is varied with the heights of the acoustic peaks 
held fixed, because Cf"^^ w -C^"^". 



VI. SIMULATIONS 



In this section, we will demonstrate our quadratic es- 
timator in Monte Carlo simulations which include the 
patchy reionization contribution to the CMB. Fully re- 
alistic simulations are outside the scope of this paper; 
we will make some simplifying assumptions, deferring a 
more complete study to future work currently in prepa- 
ration [49] . Our most critical simplifying assumption will 
be to simulate the lensing component of the B-modc as if 
it were a Gaussian field with the power spectrum of the 
true lensed B-mode. 

In fact, the true statistics of the lensed B-mode are 
non-Gaussian and this is a significant concern for our 
quadratic estimator. Our power spectrum estimator CJ'' 
is a four-point estimator in the GMB, and gravitational 
lensing generates a large four-point function [4, 40, 41], 
which could potentially be a source of bias. One promis- 
ing solution to this problem is to apply delensing algo- 
rithms [42] : if the four-point function induced by gravita- 
tional lensing is large enough to be detected statistically, 
then the level of lensing contamination can be reduced 
by partially reconstructing the lensing potential. This 
procedure can be iterated [9] to obtain an optimal re- 
construction which minimizes the level of the residual 
four-point function. Preliminary results from work in 
preparation indicate that an analysis procedure which 
combines delensing with the reionization estimator from 
this paper, by solving simultaneously for the (f> and At 
fields, is a promising approach, but this is outside the 
scope of this paper and so we use Gaussian simulations 
here. 
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A. Practical estimator 

The harmonic-space form of the estimator given in §IV: 



Tin 



N, 



V 



(60) 



^SB(r) |2 



is convenient for forecasting signal-to-noise, but too 
sfow for practical evaluation: the computational cost is 
^(^max)i which is prohibitively expensive for ^max ^ 200. 

Fortunately, Eq. (60) can be written in a fast position- 
space form. First we define spin-2 fields x ^-nd -0: 



X(n) 
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af„(2y,„(n)) (62) 



af„(2r,„(n)) (63) 



Then the estimator Tf„i is given in position-space form 
by 



ci2ny;„(n)[x(n)Xn)-fx(S)0(n)*] (64) 



In this form, the computational cost of the estimator is 
^(•^max)i using a fast isolatitude spherical transform. 

It is also convenient to have a fast position-space form 
for computing the noise power spectrum NJ"^ . If we de- 
fine correlation functions and Cj; by 



Cx" - ±E 
c." = E 



2^-1-1 
An 

21+1 



(0) (65) 



47r 



(Cf'^O^ 



££"2, ±2 



{0), (66) 



then A'J'^ is given by 



{Nrr' = ^ J_^d{cos0)dUo)[ci{o)c^{o)+c^m^m, 

(67) 

where ^2 -1-2 ^-^d doo reduced Wigner D-functions. 
This form of NJ^ is faster (e'(^^ax) vs. 0(£^ax)) than 
the harmonic- space form given previously in Eq. (61). 

Before presenting Monte Carlo results, there is one 
more ingredient. The most straightforward way to es- 
timate a T-power spectrum would be to use the power 
spectrum estimator: 



CJ 



2i 



-E 
+ 1^ 



(68) 



where Tim is the quadratic estimator. Note that this 
definition of the estimator includes subtraction of 
the noise bias NJ'^ . In turn, NJ'^ is a function of the 
CMB power spectra Cf^ and Cf^. One can ask, if 
the observed C^'s in our sky can differ (within cosmic 
variance) from the model C^'s, should we compute the 
noise bias using the model Ci& (as we have assumed when 
defining in Eq. (68)) or the observed C^'s? 

In fact, simulations show that the variance of the esti- 
mator CY is much smaller if the noise bias is recom- 
puted in each Monte Carlo realization, using the observed 

's in the realization. Intuitively, this reduction in vari- 
ance can be thought of as subtracting scatter in the noise 
power spectrum N'l'^ which arises because the CMB CV's 
will fluctuate (within cosmic variance) from realization 
to realization. Formally, we define an improved power 
spectrum estimator Ci'"^' in the following way: 



= 2m ^ ~ ^^^^ 



where the estimator NJ^ is defined by: 



rr _ {W? 



21- 



^(c,';E + jv,'^E)^(cB= + iv«fl)^ 



E* E B* B 

(24 + i)(2£2 + i) 



(70) 



Note that a fast position-space algorithm for computing 
A^J'^, analogous to Eqs. (65)-(67), is given as follows: 



nB*nB 



4^f- (Cf^-f iVf^)2"2,±2 
/(-,EoEi^2 E* E 



4,±2(^) (71) 



Nr = <Nrf j ^d{cos0)dUd)Km;{e)+Q{0)c^m 

In addition to the improvement in variance, another 
advantage of C™^' is that the noise bias subtraction does 
not depend on having precise knowledge of the fiducial 
model. The model CiS are only used in the estimator 
(via the {Ci + Ng) denominators in Eqs. (60), (61) and 
(70)) for weighting purposes. The improved estimator 

C™'"^ only depends weakly on the choice of model C^'s 
and, more importantly, is not biased if the wrong model 
is used. ^ 

We have defined C™^^ in such a way that the ex- 
pectation value (C""^') is always zero when taken over 
any ensemble of Gaussian CMB realizations (even if the 
CMB power spectra differ from the fiducial ones), but 
becomes nonzero in the presence of the non-Gaussian 
patchy reionization signal. This is particularly impor- 
tant for the reionization E-mode bump at low £, where 
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the model C^'s may have large uncertainties due to the 
number of reionization parameters that must be fit from 
the reionization bump itself. 



B. Monte Carlo pipeline 

The first step in our Monte Carlo pipeline is to simulate 
the observed CMB polarization [Q ±iU]ohs- Our simula- 
tion procedure includes several simplifying assumptions 
as we now describe. We simulate At as a Gaussian field 
with power spectrum given by Eq. (6), rather than simu- 
lating bubble formation during reionization. We simulate 
the primary E-mode Eq and response field Ei as corre- 
lated Gaussian fields, making the constant quadrupole 
assumption (from §IV) so that Ei can be treated as a 2D 
field. In a more complete treatment we would include the 
redshift dependence of Ei as described in §VB. As dis- 
cussed above, we simulate the lensed B-mode -Bicnsod as if 
it were a Gaussian field with the correct power spectrum, 
rather than including the true non-Gaussian statistics. 
We then compute the observed polarization: 

[Q±iC/]obs(S) = [Q±^C/]o(S) + [Q±^C/]lc„scd(S) 

+ {AT){n)[Q±iUUn) (72) 

Note that these simulations do not include noise. We will 
present Monte Carlo results for an ensemble of cosmic 
variance limited simulations to ^max = 2000. 

We next apply the power spectrum estimator C™'"^ 
defined in Eq. (69) to the observed polarization in each 
realization. 

Results from 1000 Monte Carlo simulations are shown 
in Fig. 12. We bin the estimated power spectrum in 
bands with A£ = 200 and we assign a bandpower that is 
the average of £{£ + l)C™'"^/(27r) over each band. In an 
ensemble of simulations with homogeneous reionization, 
it is seen that the estimated r-power spectrum has mean 
zero, and bandpower errors which are consistent with 
calculations in §V. When patchy reionization is included 
in the simulations, the estimated r-power spectrum ac- 
quires a nonzero expectation value which matches the 
fiducial power spectrum at low £, but is biased high 
at high £. 

This power spectrum bias is a known phenomenon in 
the context of lens reconstruction, where the expectation 
value of the estimator Cf'^ is linear in the power spectrum 
Cf"^ of the lens potential, but the matrix which relates 
the two is not the identity matrix [43] : an additional term 
appears when one writes out all possible contractions us- 
ing Wick's theorem. The bias can be absorbed into the 
normalization of the power spectrum estimator, either 
by computing the relevant matrix analytically, or using 
an iterative approach as suggested in [43]. However, we 
will not characterize the bias in detail in this paper; our 
objective here is simply to demonstrate that our signal- 
to-noise forecasts in §V can be achieved in simulation, 
under some simplifying assumptions that we will relax in 
future work. 
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FIG. 12: Comparison between forecasted and simulated r- 
power spectrum errors, in bandpowers with A£ = 200 and 
assuming cosmic variance limited measurements up to ^max ~ 
2000. Black/left: Forecasted power spectrum errors, us- 
ing the signal-to-noise calculations from §V. Blue/middle: 
Monte Carlo power spectrum errors using the estimator C™'"^ 
(Eq. (69)), in an ensemble of simulations without patchy 
reionization (r = 0.084, Ay = 0). Red/right: Monte Carlo 
power spectrum errors in the fiducial patchy reionization 
model (r = 0.084, Ay = 19.0). The dashed curve is the fidu- 
cial power spectrum CJ^. 



VII. FORECASTS 

We have now constructed a quadratic estimator Tim 
for fluctuations in the optical depth due to patchy reion- 
ization, shown how to compute signal-to-noise (§V), and 
tested our signal-to-noise calculation against simulations 
(§VI). In this section we will make more detailed fore- 
casts. All results in this section use the simplified esti- 
mator from §V C. 



A. Signal-to-noise 

In Fig. 13, we show the total signal-to-noise of the r 
reconstruction for a range of noise levels and beam sizes. 
A satellite experiment with /sky = 0.7 and ^fwhm ~ 4 
arcmin can obtain a 3cr detection if the instrumental 
noise is < 4 /iK-arcmin. For a ground-based experiment 
with /sky = 0.05 and ^fwhm = 1 arcmin, the thresh- 
old for a 3tT detection is Ap = 0.7 ^K-arcmin. We also 
forecast total signal-to-noise using parameters for the fu- 
ture SPTpol and EPIC experiments (Table I) and find 
(S/N)^ = 0.3 and (S/N)'^ = 28 respectively. We note 
that it may be possible to improve these signal-to-noise 
estimates using delensing or cross-correlating T£,„ with 
large-scale structure; see §VIII for more discussion. 

The preceding forecasts assume that the number of 
modes that can be measured at a given angular scale is 
proportional to /sky However, a real experiment using 
differential detectors on a small patch of sky cannot mea- 
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V 


Sfwhm 


Ap 


/sky 


SPTpol 


90 GHz 


1.7 arcmin 


4.3 /iK- arcmin 


0.015 




150 GHz 


1.0 arcmin 


4.0 /iK- arcmin 


0.015 




220 GHz 


0.7 arcmin 


12.0 /iK-arcmin 


0.015 


EPIC-2m 


30 GHz 


26.0 arcmin 


19.20 ^K-arcmin 


0.7 




45 GHz 


17.0 arcmin 


8.27 /iK-arcmin 


0.7 




70 GHz 


11.0 arcmin 


4.19 /iK-arcmin 


0.7 




100 GHz 


8.0 arcmin 


3.24 /iK-arcmin 


0.7 




150 GHz 


5.0 arcmin 


3.13 /iK-arcmin 


0.7 




220 GHz 


3.5 arcmin 


4.79 /iK-arcmin 


0.7 




340 GHz 


2.3 arcmin 


21.59 /iK-arcmin 


0.7 



power spectrum ISSJ'^ assuming Ap = 0.5 ^K-arcmin, 
dpwHM = 1 arcmin and /sky = 1- The ratio of the two 
gives the signal-to-noise per mode, which is < 0.025 for 
all multipoles, while the overall detection is ~ 15cr. 

We can ask: how do different angular scales in the 
CMB contribute to the signal-to-noise of the r recon- 
struction? To answer this question, we write the inverse 
noise power spectrum as a sum over "E multipoles" is 
or over "B multipoles" Ib'- 



NJ 



(73) 



TABLE I: Experimental parameters used when making fore- 
casts for SPTpol and EPIC in §VII. In each case, we assumed 
that the lowest and highest frequency channels were used for 
foreground subtraction, i.e. forecasts were performed using 
the sensitivity of the 150 GHz channel alone (SPTpol) [44] or 
the middle five channels (EPIC), taken from [45]. 
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(74) 



sure the reionization E-mode at ^ w 8, since the wave- 
length will be larger than the patch size. For a small- /sky 
experiment, the preceding forecasts will apply without 
modification if a cosmic variance limited measurement of 
the large-scale E-mode is available from external data. 
We find that the signal-to-noise is degraded by < 40% if 
no external large-scale E-mode measurement is available, 
or < 10% if the large-scale E-mode has been measured 
with roughly the sensitivity expected for the Planck satel- 
lite. 
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FIG. 13: Total signal-to-noise of the r reconstruction as a 
function of noise level Ap and beam size ^fwhm, assuming 
^max = 4000. 

The total signal-to-noise of the r reconstruction can 
be large enough to give a many-sigma measurement of 
the power spectrum, but the signal-to-noise per mode 
of the Ar-field is always low. In Fig. 14 (top panel) 
we compare the signal power spectrum CJ"^ to the noise 



^ (75) 

When we reconstruct r fluctuations on some scale £, the 
reconstruction is a quadratic combination of E-modes 
and B-modes whose angular scales Ie and is need not 
be the same as the angular scale I of the t mode that 
is being estimated. The quantities ffif^s) and ffiis) 
defined in Eqs. (74) and (75) can be interpreted as the 
distribution of statistical weight in £e and £b for a fixed 

e. 

In Fig. 14 (middle and lower panels), we show fPils) 
and ffiis) for three different values of £. It is seen that 
the statistical weight is a sum of two contributions. The 
first contribution appears as a "bump" with ^ 30 
and a "spike" at £b ~ £■ This three-way correlation 
between {E,B,t) can be interpreted as modulation of 
the large-scale reionization E-mode by a smaller-scale r 
fiuctuation, generating a B-mode with nearly the same 
wavelength as the r fluctuation. (Note that if £b is small, 
then £b must be « ^ because the triple (^_e,^b,-^) must 
satisfy the triangle inequality.) The second contribution 
comes from a wide range of multipoles £e and £b with a 
shape which is roughly independent of £. This contribu- 
tion can be interpreted as arising from the e""^ screening 
of the acoustic peaks, which generates B-modes across a 
wide range of angular scales. 

These two regimes can be considered separately, and 
one can construct two quadratic estimators Tim- one 
coming from a region with £e < 50 and the other one 
coming from £e > 30. The two estimators could be cross 
correlated in order to have a more robust power spectrum 
estimate, or they could be subtracted from each other, 
as a null test. We find that the detection significance ob- 
tained using this cross-correlation is w 70% (in the limit 
of low instrumental sensitivity) to w 40% (in the limit of 
high sensitivity) of the total detection significance. 
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function £[7ri |data] to the entire parameter space by com- 
puting second derivatives at maximum likelihood, and 
assuming that the likelihood function is a multivariate 
Gaussian. This should be a good approximation in di- 
rections of parameter space which correspond to well- 
constrained combinations of parameters, but in poorly- 
constrained directions the approximation may not be ac- 
curate. For example, in a parameter space which con- 
tains degeneracies, the la region of parameter space can 
have a complex shape but the Fisher forecast will always 
approximate it by an ellipse-shaped region oriented in a 
direction which is tangent to the true region. In this situ- 
ation, the Fisher matrix is still useful for identifying the 
presence of a degeneracy, and the corresponding direc- 
tion in parameter space, but the forecasted uncertainty 
in the degenerate direction is not necessarily accurate. 



FIG. 14: Top panel: Comparison between signal and noise 
per mode of the Ar-field, assuming i'max = 4000 and an ex- 
periment with Ap = 0.5 /iK-arcmin and Ofwhm ~ 1 arcmin. 
Middle panel: Distribution of statistical weig ht fFiis), de- 
fined in Eq. (75), in the B-mode multipole £b for r fluctu- 
ations reconstructed at £ = 150, 400 and 1500, showing a 
"spike" aX £b = £ on top of a broad distribution. Bottom 
panel: Distribution of statistical weight ffits) in the E-mode 
multipole £e, showing a bump at low Ie plus a broad distri- 
bution that extends to high £e. 
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B. Reionization model parameters 

What can we learn about patchy reionization by mea- 
suring the CJ"^ power spectrum? As a first look at this 
question, we will forecast parameter uncertainties in the 
five-parameter model {R^ fJiaR, r, Aj,, 6} described in §11. 

This simple semianalytic model does not contain some 
of the qualitative features seen in simulations, such as 
the dependence of the characteristic size of the ionized 
regions on the redshift [1] , or the dependence of the bub- 
ble bias on the size of the ionized bubbles [25]. However, 
we will restrict our scope to the five-parameter model to 
get a qualitative sense of the information contained in 
CJ"^. All forecasts in this section are for the EB estima- 
tor, assuming Ap = 0.5 /xK-arcmin, ^fwhm = 1 arcmin 
and /sky = 0.7. 

Our tool for forecasting will be the Fisher matrix: if 
the power spectrum CJ'^ is a function of N parameters 
{tti, . . . , TT^v}, then we compute the N-hy-N matrix as 

- —1.(2^+ 1) (c- + 7V-)2 (76) 

The rms uncertainty on parameter tt^ if the remain- 
ing [N — 1) parameters are assumed fixed is given by 
a{TTi) = {Fii)~^/^; if the remaining parameters are as- 
sumed marginalized then the uncertainty is (7(77^) = 

Some caveats should be kept in mind when interpret- 
ing Fisher matrix forecasts. The Fisher forecasting pro- 
cedure can be interpreted as extrapolating the likelihood 
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FIG. 15: Forecasted uncertainties on the total optical depth 
r and reionization width parameter Ay, assuming that the 
remaining parameters of the model {R,ainR,b} are fixed to 
fiducial values. Ellipses here and throughout § VII B are plot- 
ted at Ax^ = 1 and not 68% C.L. 

As a first step toward understanding the five- 
parameter space, consider the forecasted uncertainties on 
the parameters {r, Aj,} with the remaining model param- 
eters {b, R,ainR} fixed to fiducial values. Each of r and 
Ay can be constrained to ~ 11% if the other one is as- 
sumed fixed. The error ellipse in the Fisher approxima- 
tion is shown in Fig. 15. It is seen that the parameters r 
and Ay are nearly degenerate (the correlation coefficient 
is —0.985). The interpretation is that only one combina- 
tion of r and Ay is observable. To a good approxima- 
tion (see Fig. 16), the observable parameter is simply the 
overall amplitude A of the r-power spectrum, which we 
define by: 

^ = y^^-e(.)[l-.e(.)] (77) 

To motivate this definition of A, note that if Pax^Ax^ is 
a slowly varying function of x ^nd k in Eq. (6), then CJ"^ 
will depend on the reionization history only through an 
overall factor of A. (As a technical point, note that this 



16 



4x1 0" 



3x10" 



2x10" 



1x10" 



At = 0.084, A,= 19.0, A = 0.07 

_ _ At = 0.1 1, A,= 14,7, A = 0.Q7 
At = 0.1 1, A,= 19.0. A = 0,09 



1 000 



2000 



3000 



R = 5 Mpc. (j|„j = 0.693, R,„ = 34.1 Mpc _ 

__R = 6 Mpc, ir|„g = 0.659, R,„ = 34.1 Mpc I 
R = 6 Mpc, (7|„s=0.693, R,„ = 41.0 Mpc 




2000 



3000 



FIG. 16: Power spectrum CJ^ for different values of the opti- 
cal depth T and reionization width parameter Ay. Changing 
either r or Ay produces a change in the overall amplitude of 
the power spectrum, but if these parameters are varied jointly 
in a way which preserves the value of A (defined in Eq. (77)), 
then the amplitude remains approximately fixed. 



FIG. 18: Power spectrum CJ^ for different values of the mean 
R and width uinij of the log- normal bubble radius distribu- 
tion. Changing either R or aini? will change the location in 
I where the r-power spectrum peaks, but if the two parame- 
ters are jointly varied to preserve the value of R^s (defined in 
Eq. (78)), then the peak location is approximately unchanged. 



integral formally diverges as z ^ 0, because Xe does not 
precisely equal 1 in our fiducial model, but x(.^) — > in 
the denominator. We took Zmin = 0.1 and tested that the 
value of the integral is stable over the range [0.01,1.0].) 
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FIG. 17: Error ellipse in the Fisher approximation for the 
bubble bias b and total optical depth r, with the remaining 
parameters of the model {R, crinH, Ay} assumed fixed to fidu- 
cial values. 

Considering next the bubble bias 6, we find that b and 
T are nearly 100% correlated (Fig. 17). The bias can be 
thought of as a third parameter (together with r and A.y) 
which contributes to the overall amplitude of the power 
spectrum; however the dependence of the amplitude on b 
is relatively weak: a change of A6 = 1 produces a « 3% 
change in CJ'^. 

Now consider the parameters {_^, crinfl}. In Fig. 19, 
we show the error ellipse in the {i?, (J\nR\ plane, assum- 
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FIG. 19: Error ellipse for the characteristic bubble size R and 
dispersion ai-nR, assuming that the remaining parameters of 
the model {r, Ay, fe} are fixed to their fiducial values. 



ing that the remaining parameters {r, Ay, 6} arc fixed to 
their fiducial values. The parameters R and a\nR are 
degenerate but each can be constrained to k, 8% if the 
other is assumed fixed. This is completely analogous to 
the (r, Aj,) degeneracy described previously: there is only 
one observable combination of the two parameters. Here, 
we interpret the observable as the peak location in (., 
which we define by: 



RcS = Re 



(78) 



This definition follows [29]: Rcs is the scale where the 
square of the averaged window function {W{kR)'^) makes 
a transition from low to high values of k. The peak loca- 
tion degeneracy is illustrated directly in Fig. 18. 

In conclusion, the Fisher matrix analysis in this 
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section has shown that in the five-parameter space 
{R,ainR,T, Ay,b}, there are only two independent ob- 
servables: first, an "amphtude" observable A which is as 
a function of {r, Aj,, b} and second, a "peak location" ob- 
servable i?off which is a function of {i?,CTini?,}- The two 
observables are correlated but not degenerate: wc find a 
correlation coeflacient « —0.5 between A and RcS- 

The amplitude A of the r-power spectrum only de- 
pends on the mean ionization fraction Xe{z). The mean 
reionization history can also be constrained using the 
large-scale E-mode, but A isolates the patchy epoch 
in the sense that it only acquires contributions from 
redshifts for which reionization is inhomogcneous (0 < 
Xe{z) < 1). A comparison between parameter constraints 
from the r-power spectrum and the large-scale E-modc 
is likely to depend strongly on the parameterization of 
the reionization history which is assumed. 

Another potential source of information is the small- 
scale temperature power spectrum [46] , which is sensitive 
to patchy reionization via the kSZ effect, but requires 
careful modeling of other secondaries on small scales (e.g. 
gravitational lensing, the kSZ effect from low-redshift 
large-scale structure, thermal SZ, and the Rees-Sciama 
effect). It would be interesting to compare parameter 
constraints from the r-power spectrum and the small- 
scale temperature power spectrum, but this is outside the 
scope of a first paper. One advantage of working with the 
r reconstruction is that the greatest signal-to-noise comes 
from modes on sufficiently large scales {£ ~ 400) that 
the patchy reionization contribution dominates the low- 
redshift contribution from large-scale structure. (This 
can be seen directly by comparing the two curves in 
Fig. 2.) The r-power spectrum is also weighted toward 
high redshifts (relative to the small-scale temperature) 
because the kSZ power spectrum contains a factor (w^) 
and the velocity power spectrum is growing in time. 



VIII. DISCUSSION 

In this paper we have introduced a new technique for 
extracting the patchy reionization signal from the cosmic 
microwave background. In a fixed realization of the Ar- 
field, the two-point function of the CMB acquires terms 
which are linear in Ar, and this allows us to construct 
a quadratic estimator rf,„ for the modes of the Ar-field. 
This construction is formally similar to the lens recon- 
struction estimator (pim which estimates modes of the 
CMB lens potential (j) using the induced two-point func- 
tion of the observed CMB. 

In principle the general construction in this paper can 
use all combinations of the {T,E,B} fields, and can dis- 
tinguish transverse modes of the Ar-field at different 
redshifts using the redshift dependence of the induced 
CMB two-point function. However, we find that in prac- 
tice, the redshift dependence is so weak that only one 2D 
field (Ar)^™ can be reconstructed, and that nearly all the 
signal-to-noise comes from the EB quadratic estimator 
(with contributions from Thomson and screening effects 
which are comparable in signal-to-noise). This analysis 



also implies that it is possible to write the estimator in 
a simple form (§VC) which ignores subtleties like the 
redshift dependence. 

It is worth emphasizing that the quadratic estimator 
framework depends on having a large cross-correlation 
between the new CMB anisotropy generated by patchy 
reionization and the primary CMB. This is why our esti- 
mator is sensitive to the Thomson and screening effects 
during patchy reionization: in both of these cases, the 
new CMB anisotropy is highly correlated to the large- 
scale reionization E-mode, or to the acoustic E-mode, 
respectively. However, the quadratic estimator is only 
sensitive to the kSZ effect in temperature via its cross- 
correlation to the Doppler effect on large scales. Another 
way of saying this is that the quadratic estimator can only 
"see" the (Doppler, kSZ) cross-correlation but not the 
(kSZ, kSZ) power spectrum, since the two are ©(Ar) and 
©(Ar^) effects respectively. For this reason, we consider 
it a separate unsolved problem, not addressed in this pa- 
per, to construct a statistic which can separate the kSZ 
signal in temperature from other sources of small-scale 
power such as CMB lensing or the primary anisotropy. 

In this paper we have done a complete analysis of the 
signal-to-noise that can be achieved using a quadratic 
estimator construction, concluding with a Fisher matrix 
forecast for a simple five-parameter reionization model 
(§VII B). In this model, there are parameter degeneracies 
that leave only two well-constrained power spectrum ob- 
servables: an "amplitude" observable A which is mainly 
sensitive to the duration of patchy reionization, and a 
"peak location" observable i?off which is mainly sensitive 
to the size distribution of bubbles. 

Let us conclude by commenting on some possible con- 
tinuations of this work which are outside the scope of a 
first paper. 

One such extension is to increase the realism of the 
reionization model, and to study prospects for constrain- 
ing a more complex parameter set. It is unclear whether 
the "two-observable" picture that we obtained in our 
Fisher matrix analysis of the five-parameter model will 
be qualitatively different in a more complete model. 

Another natural continuation is to study prospects for 
cross-correlation with other flavors of cosmological data 
which probe the epoch of inhomogeneous reionization. 
A nice feature of the quadratic estimator framework (as 
compared to other statistics such as the CMB power spec- 
trum) is that the patchy reionization signal is extracted 
in the form of a map rgm which can be cross-correlated to 
other datasets. This cross-correlation technique is part 
of future work, and can be used either to boost signal-to- 
noise for a first detection (e.g. [10, 11] in the context of 
CMB lens reconstruction) or to make the measurement 
more robust by reducing the range of systematic effects 
which may contaminate the auto power spectrum of the 
reconstruction. 

Finally, a continuation of this work which is cur- 
rently in preparation is the implementation of the es- 
timator in more realistic simulations which include the 
non-Gaussian statistics of the lensed B-mode. The simu- 
lations in this paper (§VI) show that our signal-to-noise 
forecasts can be achieved in toy Monte Carlo simulations 



18 



in which the Icnscd B-mode has the correct power spec- 
trum but the statistics are treated as Gaussian. Including 
non-Gaussianity is an important step: it can potentially 
bias the estimator T£,„ (which reconstructs the modes 
of the Ar-field by using the induced non-Gaussian sig- 
nature) but also presents the opportunity for improving 
signal-to-noise by combining t reconstruction with dc- 
lensing. 
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APPENDIX A: QUADRATIC ESTIMATOR 
FORMALISM 



where: 



1 



2(2^+1) 



E 



XYX'Y'iit2 



pX'y'(T)*.p,-ixYF' 



(A6) 



The purpose of this appendix is to develop the 
quadratic estimator construction that was used in 
§IV, §V and to estabhsh some key properties: the form 
of the quadratic estimator (Eqs. (36), (48)), the noise co- 
variance (Eqs. (37), (49)), and the total signal-to-noise 
(Eq. (52)). In §A 1 we will develop the "single- field" 
version that was used assuming the constant quadrupole 
approximation in §IV. In § A 2 we will generalize to the 
multifield case that was used in §V to correctly treat the 
redshift dependence. 



where Ci is defined in Eq. (50). Now that we have 
found the optimal estimator, we consider its signal-to- 
noise. Let us carefully distinguish between what we mean 
by "signal" and what we mean by "noise" , in the context 
of T reconstruction. 

We write the two-point function of Tim as the sum of 
two terms: 



1. Quadratic estimator formalism: single field 

Let us first make the approximation that patchy reion- 
ization is narrow, so that the anisotropic optical depth 
can be described by a single 2D field (Ar)^„i. In §V, 
we showed that to lowest order in At, patchy reion- 
ization produces correlations between CMB multipoles 
flfm 7^ oii'm' which are of the form: 



mi m2 m 



^^^(Al) 

where X,Y G {T,E,B}, and the couplings F^^^/-* are 
given in Eqs. (42)-(46). 

Given the two-point function in Eq. (Al), we construct 
a quadratic estimator Tim for each mode of the Ar-field, 
by solving for the weights Wi^^^„^^,,„, in 

XY iimil2m2 

which minimize the variance 

Var(T/^) — - ^ ^ ^hn^imii2m2 



XYX'Y' £imi£2ni2 
^ " e'm'eimil2m2^ei 



C^.'^'CZ, (A3) 



subject to the normalization constraint (r^m) = (Ar) 



KTlirnTl'm') 
{'''im'^^'m') noise 



{^em'''e'm') noise + ('''frn'^^'m' )signal 

crs (A7) 



The first term (r? T£, 



m' / noise 



is obtained by summing 



Gaussian contractions between CMB multipoles af^^ and 
is what we mean by "noise" : it is the power spectrum 
of the reconstruction due solely to statistical fluctuations 
in the CMB, in the absence of any patchy reionization 
signal. The second term (T/m^f'"i')signai arises from the 
expectation value (jim) = (Ar)£,„ and is what we mean 
by "signal" : excess power in the reconstruction due to the 
presence of r fluctuations. We therefore interpret NJ'^ as 
the reconstruction noise power spectrum and CJ'^ as the 
signal power spectrum. 

In a patchy reionization model with signal power spec- 
trum CJ'^, the total signal-to-noise of the reconstruction 
(summed over all modes) is 



S/N 



/sky 

2 



^(2£+l) 



CI 



NJ 
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A short calculation shows that the minimum-variance es- 
timator is given by: 
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><{C-')ff afX{C-^)jr aJX^_. (A5) 



This is the "number of sigmas" for an overall detection 
of patchy reionization, via excess power in the t recon- 
struction. 

Our construction of Ti„i is formally identical to the 
quadratic estimator for the CMB lens potential 4>. 
More precisely, to lowest order in (p, the two-point func- 
tion induced by lensing is of the form of Eq. (Al), with 

the Ar-field replaced by 0, and the object ^f^^P re- 
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placed by 
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TE ^ i^tit. 



-iTTpO 
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We propose two possibilities for constructing a 
quadratic estimator in the multifield case. The first pos- 
sibility (which corresponds to the "simple estimator" in 
§V C) is simply to use the single-field quadratic estimator 
corresponding to a fixed redshift bin /x: 



'em 



n: 
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The only new ingredient in the multifield case is the non- 
(A9) trivial redshift response: the expectation value is given 
by 



where the F symbol is defined by: 

= [-^1(^1 + 1) + ^2(^2 + 1) + 4(4 + 1)] 



^ / (2£i + l)(2^2 + l)(2£3 + l) / h I2 4 U^^Q. 

167r \-s s ' 



With these replacements, it can be checked that the es- 
timator we have derived in Eq. (A5) agrees with the well- 
known minimum-variance quadratic estimator for 
lens reconstruction. In the single-field case, the quadratic 
estimators for lens reconstruction and patchy reionization 
can simply be regarded as arising from a different set of 
couphngs r^^^. (See also [47] for another example of the 
quadratic estimator formalism: reconstructing a spatially 
fluctuating rotation a(n) in the CMB polarization.) 



where 



rT(/j) 



2(2£ 
xT 



' XYX'Y'lil2 
X'Y'{T„)*,f^-l.YY' 



Y.XY{t^) ,^-l^XX' 
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(Note that Rf\z>') = 1.) 

It follows that the signal power spectrum of the "sim- 
ple" quadratic estimator is given by 

(?i™*4™'>.ignal = Y.^R['\z'')fCl-^-6u'5mra' (A16) 



2. Quadratic estimator formalism: multiple fields 

The quadratic estimator formalism in the preceding 
subsection can be applied whenever the CMB two-point 
function {af^(i{,^i) is proportional to an auxiliary field 
which is not directly observable and must be estimated, 
such as Ar or <j). (The form of the two-point function 
given in Eq. (Al) is the most general form allowed by 
global rotation invariance.) 

However, in §V we will need greater generality: we are 
also interested in the case where there are N auxiliary 
fields At", where a = 1, . . . , N runs over redshift bins. 

We consider a CMB two-point function of the general 
form 



i«.,™,«.2™2) - 2.^ w km, 

(All) 

for some set of couplings ^f^J^l" ' • This form of the two- 
point function keeps track of the (weak) dependence on 
the redshift bin. 



and the total signal-to-noise is given by: 
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/sky 

2 



^(2£+l) 



N 



r(M) 



1/2 



This construction has the advantage of preserving the 
simplicity of the single-field estimator in the multifield 
case, but the disadvantage that it does not necessarily 
extract all the signal-to-noise. 

The rest of this appendix is devoted to a more com- 
plicated eigenmode construction (this corresponds to the 
principal component analysis in § V B) which is math- 
ematically guaranteed to contain all the signal-to-noise. 
For an arbitrary set of couplings T^^^/"-*, the total (S/N) 
of the simple estimator (Eq. (A17)) will be less than 
the total (S/N) of the eigenmode estimator (Eq. (A26)). 
However, in our fiducial reionization model, the calcula- 
tions in §V C show that the signal-to-noise is almost the 
same in the two cases. Therefore, in practice the eigen- 
mode construction will serve as a proof that the simple 
construction does not lose information. 
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For each £ we define an N-hy-N matrix 



f: 



af) 



2(2^+1) 
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For each £, we define weights by solving the 

eigenmode equation: 



/3 



(A19) 



where the eigenvalues A^'^ satisfy A^^^ > > ■ ■ ■ > 

AN) 

We then define an eigenmode estimator by: 
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E E 
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We also give the explicit form of the Ti-ficld, which was 
left as a functional derivative in Eq. (24) and not evalu- 
ated explicitly. 

We will work in the synchronous gauge where h and 
T] denote the scalar degrees of freedom of the metric. 
Overdots denote derivatives with respect to conformal 
time. We denote the visibility function by g{x) = Te~'^ . 
Further notation in this appendix follows [48]: is 
the baryon velocity, I^ti and Ap.f denote the Lcgcn- 
dre expansion of the local temperature and polariza- 
tion, respectively, 11 = At2 + ^po + ^P2 is the source 
for polarization generated by Thomson scattering, and 
a= (/i-H6?7)/(2fc2). 

Let us begin by calculating Cf°^'^. The hne-of-sight 
integral for polarization (Eq. (15)) can be rewritten in 
harmonic space [48] as: 



xy'dxff(x)^^n(x,k) (Bl) 



af3 



Because the perturbations are linear, n(x, k) depends 
(A21) linearly on the initial curvature perturbation C(k) and so 
we can rewrite Eq. (Bl) in the form: 



The mean response of the eigenmode estimator is given 
by: 

(^™) = (A22) 

a 

Rf\z") = Nl^'^'^Y.K^'^fi^^) (A23) 

a 

and a short calculation shows that the signal and noise 
power spectra are given by: 

(^lm*^''m') noise = N^^^^^ SijSu'^mm' (A24) 
(^2*41')-gnal = X^;^N;^^'^S,j6u,Sram' (A25) 

Because estimators i ^ j arc uncorrelated, it is straight- 
forward to calculate total signal-to-noise: we simply sum 
over eigenmodes. 
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The eigenmode construction is convenient for study- 
ing signal-to-noise properties in the multifield case: one 
can compute the total S/N achievable in principle with 
quadratic estimator, or count the number of eigenmodes 
with significant S/N in the sum in Eq. (A26). 



= Am' J d^k Af (fc)C(k)F;„(k) (B2) 

This representation encodes all the evolution in a single 
transfer function Af and will be convenient for comput- 
ing power spectra below. 

The response field Ei can be treated analogously. We 
rewrite the linc-of-sight integral for Ei (Eq. (18)) in har- 
monic space: 



{Am') 
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and define a transfer function Af so that: 



(B3) 



<^t^'^ = 4^*' / ci^k Af (^)(fc)C(k)r/„(k) (B4) 



The transfer functions Af (fc), Af ^''^(fc) can be com- 
puted using CAMB, and Cf°^^ is given by: 



APPENDIX B: POWER SPECTRUM 
CALCULATIONS 



EoE,(z) _ (4^)2 y k^dkAf{k)Af''^'\k)P^{k) (B5) 



In this appendix, we will describe our procedure for cal- 
culating the power spectra Cf°^^ , where X,Y € {T, E}. 



This concludes our calculation oiCf°^^\ let us now dis- 
cuss the temperature field. 
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by 



The line-of-sight integral for temperature [48] is given where: 



(W) I d^kY;„,(k) J dxjiikx) 

x|^e-"(^)Qo(x,k)+g(x)Qi(x,k) 

+.9(x)Q2(x,k)+.g(x)g3(x,k) I (B6) 



where 
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Qi 
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(Note that the quantity in parentheses in Eq. (B6) is 
what we called St previously in Eq. (22).) 

To get Ti from T, we take the functional derivative 
(see Eq. (24)): 
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There is a subtlety here: when we take the functional 
derivative, we only differentiate the "explicit" depen- 
dence of r(n) on t(x) through the factors of e^'^, p, 
g and g in Eq. (B6). We treat the source terms Qi as 
independent of r, but this is actually an approximation 
since r does appear in the evolution equation for these 
quantities. This approximation is equivalent to dropping 
some double-scattering terms which arc suppressed by 
one power of r; we checked directly that the approxima- 
tion is good to a few percent. (The same approximation 
was made previously in polarization.) 

To compute power spectra, we define transfer functions 

AT and A?^"'-' so that 



Ti{z) 



As in the polarization case, we compute the transfer func- 
tions using CAMB. The power spectra which contain 
temperature are then given by: 
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